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A simple closed-form expression for the number of days between the winter sol- 
stice and the latest sunrise is derived. Formulas for the summer solstice and the 
sunset are derived as well. These approximate formulas make it easy to see the 
functional dependence on the latitude of the observer. An approximate expres- 
■ sion for the difference in time of day between the latest sunrise and the sunrise 

' on the solstice is also derived. The formulas are not valid in the Arctic, Antarctic 

^ . or tropics. 

O . 

^ ■ Introduction 

o 

^ . In this paper, we present a simple closed-form expression for the number of days between 

^! the winter solstice and the latest sunrise. This formula is the result of a lengthy calculation 
O I which simplifies greatly in the end. A few approximations are made in the derivation of the 
' result, but it agrees well with actual data for sunrise times. 

• Let us begin by recalling the reasons why the latest sunrise does not occur on the shortest 

day of the year. There are two effects that contribute to this phenomenon. The dominant 
effect is that due to the inclination of the Earth's axis of rotation. The other contribution is 
^ I due to the eccentricity of the Earth's orbit. 

We first discuss the effect due to the inclination of the Earth's axis of rotation, that is, 
we will begin by taking the Earth's orbit to be a perfect circle. It is useful to work with the 
concept of the celestial sphere. In the course of a year the sun's path is a great circle on the 
celestial sphere. This great circle is inclined at an angle of 23.44° relative to the celestial 
equator. The sun moves along this great circle, the ecliptic, at a constant rate, because we 
are assuming that the Earth's orbit is a perfect circle. Thus, the right ascension of the sun 
does not change at a constant rate. (The right ascension is the east/west coordinate on the 
celestial sphere that plays the same role as the longitude coordinate on the surface of the 
Earth.) At a solstice, the right ascension of the sun changes the fastest because the ecliptic 
is closest to one of the poles there. At an equinox, the right ascension changes the slowest 
because the ecliptic intersects the celestial equator at an angle. 

Thus, if a picture is taken of the sky at 12:00 noon every day for a year, we see not only 
the sun's familiar north/south motion, but also an east/west motion. The reason is that 24 
hours is the length of an average solar day, but we have seen that during the year the sun's 
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right ascension changes at a non-uniform rate. The images taken of the sun would trace out 
a certain figure eight, called the analemma [1]. When drawing the numeral eight with a pen, 
we see that the pen moves in the same direction at the top and bottom, and in the direction 
opposite to this during both passages through the center. This corresponds to what we said 
above about the solstices versus the equinoxes. 

In particular, we are interested in the fact that the sun is moving in the east/west direction 
in this sequence of photos, at the winter solstice. In the northern hemisphere, the day (in 
the sense of the amount of time that the sun is visible above the horizon) is shortest on the 
winter solstice, but because of the east/west motion of the sun, the times of both the sunrise 
and the sunset are still changing at this time of year (by the same amount, because their 
difference is minimal at this time). Thus, the latest sunrise is not on the shortest day of the 
year. 

Next, we include the eccentricity of the Earth's orbit into the discussion. The path of 
the sun on the celestial sphere is not changed. What is different is that the sun is no longer 
moving along this path at a constant rate. Thus, the rate of change of the right ascension 
of the sun is changed. The result is that the analemma discussed above becomes lop-sided, 
and this modifies the amount of east/west motion at the solstices. 

Statement of result 

Let n denote the number of days between the winter solstice and the latest sunrise. Then 
at latitude A, assumed to be not too close to the equator or the poles as explained below, n 
is given approximately by 

/ 27r \ ^cos(2q;)+cos(2A) 

n = \ o — —— cos a - — -1= , 1 

\ Ny ) v^^Hana sinA ^ ' 

where 

Ny = 365.2564 (2) 
is the number of ephemeris days per sidereal year [2], and 

a = 23.44° (3) 

is the angle between Earth's rotation axis and the normal to the ecliptic, and 



iVy(l -e cos 7/^0)2 ^ ^ 

is the change in the angular position of the Earth about the Sun per ephemeris day at the 
winter solstice [3]. In this expression, 

e = .01672 (5) 
is the eccentricity of the Earth's orbit, and 

i^o = l-927r (6) 
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is the eccentric anomaly at the winter solstice. It should be noted that this value for ipo 
is approximate; its exact value does not affect the result for n significantly. The value 
of i/jo can be determined from the date of the perihelion if one assumes that the Earth's 
orbit is an ellipse. In actuality there are perturbations to the orbit, and the dates of the 
perihelions vary quite a bit. According to the Astronomical Applications Department of the 
U.S. Naval Observatory (http://aa.usno.navy.mil/data/docs/EarthSeasons.html), the dates 
of the perihelions in the years from 1993 to 1998 were Jan. 4, Jan. 2, Jan. 4, Jan. 4, Jan. 2 
and Jan. 4. 

We see that the expression for n given in Eq. (1) becomes larger as we approach the 
tropics and that it goes to zero as we approach the Arctic or Antarctic Circle. This is 
because the cos(2q;) + cos(2A) in the numerator goes to zero as we approach the Arctic or 
Antarctic Circle. The formula is not valid for values of A that are north of the Arctic Circle 
or south of the Antarctic Circle because the sun does not rise and set on a daily basis in 
those regions. The other dependence on the latitude A is the sin A in the denominator. This 
formula is not intended for use in the tropics, as explained later in this paper, so there is no 
problem with a divergence. 

The conventions that we use regarding the sign of A are that positive values are used in 
the northern hemisphere and negative values are used in the southern hemisphere. In the 
northern hemisphere, the values for n are positive, indicating that the latest sunrise is after 
the winter solstice. In the southern hemisphere, the values for n are negative, indicating 
that the extremal (earliest) sunrise is before the December solstice. The term "December 
solstice" can be used to avoid confusion in these considerations. 

We have also derived formulas for the June solstice and the sunset. The results may be 
summarized as follows. For sunsets the value of n is always minus the corresponding sunrise 
n at the same location and for the same solstice. For the June solstice sunrise, the formula 
for n is minus the formula in Eq. (1), with the appropriate 5 value used. This 5 value can be 
found from Eq. (4) using the value of the eccentric anomaly for the June solstice, which is 
approximately tt less than the value at the December solstice. These statements regarding 
the different formulas for n are made under the assumption that the angle of inclination a 
is not greater than 90°. This is true for most planets in our solar system. 

Examples 

Table I compares the approximation given in Eq. (1) with tabulated values from the 
Astronomical Applications Department of the U.S. Naval Observatory 

(http://aa.usno.navy.mil/data/docs/RS_OneYear.html) for the solstice of December 2001. 
In tables where times are rounded to the nearest minute, there will be a range of dates that 
have the latest sunrise time. The corresponding range of n values is given in the column 
labeled "actual value of n." We do not expect our approximations regarding atmospheric 
refraction and the size of the solar disk to be good near the Arctic Circle. Near the tropics, 
n is so large that the Taylor expansion to second order is not sufficient. 
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Table I 



location (city, state) 



latitude (degrees) 



approximate n 



actual value of n 



Fairbanks, Alaska 64.82 
Anchorage, Alaska 61.22 
San Francisco, California 37.77 
Key West, Florida 24.55 



2.4 
4.5 
16.2 
28.6 



2 - 5 
1 - 8 

14 - 16 
16 - 31 



Another comparison of the approximation given in Eq. (1) with actual sunrises is given 
in Fig. 1. The comparison is made for a range of latitudes, from the Tropic of Cancer to the 
Arctic Circle. 

Derivation of Eq. (1) 

Since the derivation of Eq. (1) is a lengthy calculation, we give only a summary here. The 
analysis begins by describing the motion of the Earth in a heliocentric ecliptic coordinate 
system, as described in Ref. [3]. We use the angle 9 to denote the angular position of the 
Earth relative to the Sun. This angle also gives the position of the sun in a geocentric ecliptic 
coordinate system. We define 9 to be zero at the vernal equinox, as is standard. Next, we 
transform to a topocentric equatorial coordinate system, and then to a horizontal coordinate 
system. This final coordinate system is the one an observer at a given point on the Earth 
would use, and an angle h is used to measure the elevation of objects above the horizon, 
with /i = 0° at the horizon and h = 90° at zenith. The result for the sun is 



This equation gives sin h at any instant in time. The angle is an angle that measures the 
rotation of the Earth about its own axis. The time dependence of sin h is due to the time 
dependences of 9 and 0. The angle evolves in time at a constant rate, for our purposes, 
increasing by 27r every sidereal day. The angle 9 increases by 27r every sidereal year, but not 
at a constant rate, due to the eccentricity of the Earth's orbit. 

If we ignore refraction in the Earth's atmosphere and the finite size of the solar disk [4], 
we could find the times of sunrise and sunset by setting sin h equal to zero and solving the 
transcendental equation that results for the time t. The errors made in this way should be 
more or less constant over a small part of a year, and we hope that they do not affect our 
search for the latest sunrise too much. The next approximation is to use the same value for 
9 throughout any given day. This makes it possible to solve for the times of sunrise and 
sunset using elementary functions. Again, it should be the case that these errors do not 
vary significantly over a small part of the year and thus do not affect the day with the latest 
sunrise excessively. Now that we have a formula for the time of the sunrise on any given 
day, we make a second-order Taylor series expansion about the winter solstice. This is where 
the calculation gets complicated, but with enough work the coefficients can be simplified. 
Finally, the extremum of a quadratic function can be found easily, and the result is Eq. (1). 



sin h — cos A cos 9 cos (f) + cos a cos A sin ^ sin + sin a sin A sin ^ . 



(7) 
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More details about this calculation are given in Appendix A. An approximate expression for 
the difference in time of day between the latest sunrise and the sunrise on the solstice is also 

derived there. 

We note that Eq. (1) does not behave well when a goes to zero. Apart from situations 
where a is very small, this result should be applicable to other planets. As a goes to zero, the 
quadratic term in the Taylor expansion discussed above vanishes, so a higher-order expansion 
would be necessary. 

At the beginning of this paper, wc stated that the results were intended for use at latitudes 
that were not too close to the equator or the poles. We can now clarify these limits. Using 
terminology appropriate for the northern hemisphere, we can say that the formula should 
not be used for latitudes north of the Arctic Circle because the sun does not rise and set on 
a daily basis there throughout the year. Such daily risings and settings were assumed in the 
calculation. From the point of view of mathematics, the formula gives imaginary numbers 
north of the Arctic Circle. The Tropic of Cancer is more of a soft limit. The value for n is so 
large at the Tropic of Cancer, that the second-order Taylor expansion is no longer sufficient, 
and the results become less accurate the further south one goes (see Fig. 1). 

Classroom activities and projects 

The materials for our first activity arc a globe and a measuring tape made of cloth or 
paper. The idea is to get a better understanding of the fact that the right ascension of 
the sun changes at different rates at different times of the year. The lines of longitude and 
latitude on the globe form the same type of coordinate system as the one used on the celestial 
sphere. The purpose of the globe in this demonstration is to provide a sphere with such a 
coordinate system, and we will think of this sphere as representing the celestial sphere. Place 
the beginning of the measuring tape at the point on the globe that has longitude and latitude 
both equal to zero. Then wrap the tape around the globe so that it traces out a great circle 
that is inclined relative to the equator. An angle of inchnation of about 23 degrees would 
correspond to the inclination of the Earth's axis as discussed in the Introduction, and the 
measuring tape would then represent the ecliptic. Now choose a sequence of evenly spaced 
points along the measuring tape, and find the right ascension of each point (the longitude on 
the globe). Plotting these values on graph paper will show that the right ascension changes 
the slowest when the ecliptic crosses the equator (at the equinoxes) and fastest when the 
ecliptic is furthest from the equator (at the solstices). The effect is easier to see if a larger 
angle of inclination is used. 

An extension of these ideas can be used to create a plot of the analemma. Make a table 
of the right ascension and declination (longitude and latitude) of 36 evenly spaced points on 
the ecliptic, starting at the beginning of the measuring tape and going west around globe. 
The dechnation column of the table will also be called the Y column. Next, we will create 
a third column of the table, which we will call the X column. The values for this column 
are obtained by taking the values in the right ascension column and subtracting an integer 
multiple of 10 degrees. The integer is for the first row, 1 for the second row, 2 for the 
third row, and so on. The result is that the X and Y columns give the position of the sun 
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relative to its average position at 12:00 noon (see the discussion in the Introduction). Now 
create a plot with X and Y axes. The data points should trace out a figiu^e eight. The real 
analemma is actually lop-sided, due to the eccentricity of the Earth's orbit. Some globes 
have a picture of the analemma, often placed somewhere in the Pacific Ocean where there 
is space for it. 

For our next activity, work out the value of n at your latitude using Eq. (1). Compare 
the results with values in published tables, such as those found in tide books. You can 
also use one of the the web pages of the Astronomical Applications Department of the 
U.S. Naval Observatory, http://aa.usno.navy.mil/data/docs/RS_OneYear.html , to generate 
a table of sunrise and sunset times for your location. Next, imagine hving on Mars at the 
same latitude. What would the value of n be there? It will be necessary to look up various 
orbital parameters for Mars in order to calculate this, and the length of the Martian day will 
be needed. 

The materials for our last activity are a watch and a piece of chalk. In the Introduction 
we mentioned a sequence of photographs taken of the sky every day at 12:00 noon. The 
purpose of these photographs was to record the location of the sun at the same time on 
different days. We can carry out similar measurements by recording the location of the 
shadow of the top of a flag pole or other convenient object at the same time on a sequence of 
days. If class is not in session at 12:00 noon, another time of day can be used. The greatest 
motion of the shadow from one day to the next will be seen in March and September, near 
the equinoxes. At these times, the component of the shadow's movement in the North/South 
direction will be about 0.4 degrees per day. Shadows cast by the sun are a bit fuzzy because 
of the nonzero size of the solar disk (about half a degree). Therefore, the motion can be seen 
more clearly by waiting for several days to pass. The East/West motion in which we are 
interested in this paper is large near both the solstices and the equinoxes. Its maximal value 
is on the order of 0.1 degrees per day. At the equinoxes the motion will be in the direction 
opposite to the motion at the solstices, as discussed in the Introduction. If the shadow falls 
on an area where there are tiles, an alternative to using chalk marks on the ground woTild be 
to use the coordinate system provided by the tiles to record the position of the shadow. A 
larger project would be to record the location of the shadow once a week or so, for an entire 
semester or even a year. 
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Appendix A 

This appendix gives a more detailed description of the derivation of Eq. (1). The angle 
in Eq. (7) measures the rotation of the Earth about its own axis. The value of as a 
function of time is 

Ht) = (t^^) t , (8) 



T 

sidereal ' 

where Tgidereai ~ 86164.09 seconds is the time it takes for the Earth to rotate once about its 
axis relative to the fixed stars. For our calculation, it is useful also to define 

m-{7f^^]t. (9) 

\ ephemeris / 

with Tgphemeris — 86400 scconds = 24 hours, and the difference 

A(t) = - . (10) 
Substituting + A in for (p in Eq. (7) and using some trigonometric identities gives 

sin h — Acos(0 + Q cos A + sin a sin A sin ^ , (11) 
where A and ( are defined to be 

A = ^ (cos a sin ^ sin A + cos 6 cos A)^ + (cos a sin 6 cos A — cos 6 sin A)^ , (12) 

and 

1 /cos a sin ^ cos A — cos sin A \ , , 

C = -tan-M . n ■ ^ n r • 13 

y cos a sm (J sm A + cos 6 cos A J 

The argument of the tan^^ function on the right-hand side of the equation for ( stays 
small throughout the course of the year, so there is never a question about which branch 
of the tan~^ function to use. One way to see this is to first imagine replacing cos a with 1, 
whereupon the numerator becomes sin(^— A) and the denominator becomes cos(^— A). Both 
9 and A increase monotonically from to 27r over the course of the year, so the numerator 
stays close to and the denominator stays close to 1. In actuality, cos a is not 1, but the 
functions still stay close to and 1, respectively. 

Over the course of the year, the values of modulo 27r at sunrise and sunset vary greatly, 
because cj) measures the rotation of the Earth relative to the fixed stars. In contrast, the 
values of (f) modulo 27r at sunrise and sunset do not vary by such large amounts, because cf) 
was defined to be an angle that increases by 2tt every 24 hours (mean solar day in the year 
1900). This is the reason for transforming from to 0. If we now take the values of 9 and 
A to be constant on any given day, we can solve for the times of day of sunrise and sunset. 
Expressed in terms of (p (modulo 24 hours, so is in the range from — tt to tt) these are 

7 ^ , _i / — sin a tan A sin , ^, 

= -C ± cos M J j , (14) 
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where the upper sign is for the sunset and the lower sign is for the sunrise. Note that at 
the vernal equinox, and the time t are zero, and the sun is crossing the meridian at the 
location of the observer. Thus, (p is zero near noon, and sunrise times have a negative value 
for (J), while sunset times have a positive value for 0. 

The lower sign in Eq. (14) gives an explicit formula for the time of day of the sunrise 
on each day of the year (expressed in terms of a value). We are interested in finding the 
extremum near the winter solstice. Wc use an integer k to index the days relative to the 
winter solstice [A; = on the winter solstice. A; = 1 on the next day, etc.]. Then at some time 
on day k the value of 9 is (to order k"^) 



where Ci is d9/dt evaluated at the winter solstice, times 24 hours, and C2 is SO/dt^ evaluated 
at the winter solstice, times the square of 24 hours. A similar expansion can be written out 
for A(/c), but it does not contain a k"^ term because A is a linear function of time. When 
these expressions are substituted into the formula for the time of day of the sunrise, the 
result is a function of k that is accurate to order k"^. By calculating the first and second 
derivatives of this function with respect to fc, we can calculate a quadratic approximation 
for the function, and from this find the value of k at the extremum. We use n to denote this 
extremal value of k, and the result is given in Eq. (1). There is actually a small correction 
term involving the C2 coefficient mentioned above, but this term is smaller than the result 
stated in Eq. (1) by a factor of e, so it can be neglected. 

The quadratic approximation can also be used to calculate another quantity of interest. 
We let denote the difference in time of day between the latest sunrise and the sunrise 
on the winter solstice. This time difference can be calculated by substituting the value of n 
back into the quadratic approximation for the time of day of the sunrise. The result can be 
expressed as 



Evaluating the constants, we get that At is approximately 15 seconds times n. It is interesting 
that the constant of proportionality is independent of the latitude. Sample numbers for San 
Francisco are At = 245 seconds and n — 16, with a ratio approximately equal to 15 seconds. 
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0{k)^l7r + C^k + lc2e, 



(15) 




(16) 
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The result presented in this paper is compared with actual sunrises. The curve is a plot 
of n (the number of days between the winter solstice and the latest sunrise) as a function of 
the latitude A, as given by Eq. (1). The dots are the actual values of n, plotted for latitudes 
that are integer numbers of degrees. Gray squares are also plotted for latitudes that are 
integer numbers of degrees. A darker gray square at a value of n and A means that the time 
of the sunrise on the n-th day after the solstice is within 5 seconds of the latest sunrise at 
that latitude. The lighter gray square is for a value of 20 seconds. The "actual" times for 
sunrises were computed as discussed in Ref. [4], using the standard definition of sunrise to 
be the time when the center of the sun is 50 arc minutes below the horizon. The 50 arc 
minutes take into account the size of the solar disk and the effects of atmospheric refraction. 
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